Stars point is working a little like a navigation point made with three landmarks, but in that case bearings are changed by position lines coming from stars. To have an accurate stars point it's needfull to choose stars the most rejected (the best way is with 120 degrees between each star), and do sextant's watch a little before or after the sun's falling to be able to still notice horizon line. It's important to avoid a choice on stars to close to the horizon too because of light refraction (see Snell's law). It's difficult when you are using sextant to”go down the star on horizon's line” (It's easy to confound a star with an another close to its), That's why it's easier to use sextant with head upside down (scale above), On few words “rise horyzon's line on the choosen star” choosen. The most of time it's only possible to use stars close to horyzon's point where the sun disapears, Beyond this it's hard to see clearly horyzon's line to read an accurate angle.
To write this little program i have often used an astrometry master (introduction aux coordonnées célestes, astrométrie, Jleborgne) from Toulouse university, easy to find on the web, and others things grabbed on the lane.
1) Coordonate system, definitions, motions involved
a) Celestial coordinates system
It's working with equatorial coordinates; we will locate stars with their positions against celestial equator which is finally just equator projection on celestial's sphere. In fact it's almost the same method than how to place something on earth; to locate a star we will have:
–Right ascension the most of time symbolysed by alpha letter, it's the angle read on celestial equator between vernal point and the perpendicular line with celestial equator cuting the watched object. the original point is vernal point and the angle is read from west to east, This angle is the most of time written in hours but is easily convertible in degrees (1 hour=15 degrees, 1 degree=60 minutes, 1 minute=60 secondes), in hours it's changing from 0 to 24 hours (or 0 to 360 degrees ).
– declination angle often symbolized by delta letter, angle read between the object seen and celestial equator. It can rise until 90 degrees (-90 at south pole 90 at north pole)and it can be at least 0 degrees (equator). declination angle and latitude have a lot of common points, It's positive on north hemisphere and negative on south hemisphere.
–Vernal point most of time written with gamma letter or just "g", is more difficult to understand. on celestial sphere ecliptic plan (plan with both earth and sun) and celestial equator cross their way in two points , one cross when the sun is going from north to south hemisphere, and the another one when sun is going from south to north hemisphere, so when it's spring equinox (between the 19 and the 21 March).we will see that vernal point is not motionless. equinoxes moves about 50” per year along ecliptic from east to west. vernal point's time angle (or sidereal time) is read from greenwich meridien and is written ahso or ahgy. vernal point is almost moving with constant speed, this speed is about 15.04166 degrees . we will see that to perform our little software we will need more accuracy.
For celestial coordinates J2000.0 coordinates are the origin, in fact measurements read the first of january 2000 at noon in gregorian calendar. If we write a star's name in wikipedia for exemple, the most of time numbers noticed on the right of the page are those read at J2000.0, To have the good coordinates depending of measurement's moment, it's essential to add corrections.
b) Julian calendar, Gregorian calendar
Important to understand for our program; gregorian calendar is a little like a better julian calendar's version, the origin is the same, monthes and weeks system too. The difference is about bissextial year's choices. For julian calendar, there is a bissextial year each 4 years (2008, 2012, 2016…) and for gregorian calendar, more accurate there is a bissextial year if you find other conditions more. In our case we will only use julian calendar. In calendar ther is 12 monthes filled up alternatively by 31 (january, March, mai, july, august, october, December) or 30 days (april, June, September, november). There is one exception for february with 28 days for one year “normal” and 29 days for a bissextial year, so we have 365 days for a classical year and 366 days for a bissextial year, in fact an average annual number of 365.25 days with the Julian calendar. However this duration's time doesn't fit perfectly with a real year's time duration…So we will have to cope with this mess to find the best way to make a good program ! a funny work ! or not...….
c) What does it mean “a year” ?
A year in the Julian calendar has 365.25 days. But there are several options to decide what a year is:
–sideral year : earth revolution relative to a fixed point. It takes about 365 days, 6 hours, 9 minutes and 10 seconds.
–tropical year : time between two successive passages of the sun in the direction of the vernal point. This period will be very important for our program, including the calculation of the vernal point's hour angle. It takes about 365 days, 5 hours, 48 minutes and 45.26 seconds at J2000.0. For our century it decreases about 5.23 secondes each 100 years. There is a formula useful to calculate the average tropical year duration easy to find on the web, It's :
Average tropical year =365.24219052-61.56*10^(-6)*T-0.068*10^(-6)*T^2 263*10*(-9)*T^3 3.2*10^(-9)*T^4. In 2000 it's equal to 365.2421898 days and in 2100 it's equal to 365.2421291582 days. T is time in Julian centuries since J2000.0.
–anomalistic year : time between two successive passages of the Earth at perihelion. It's about 365 days, 6 hours, 13 minutes and 53 seconds. Perihelion is the moment where earth is the closest to the sun (Kepler's laws), so it's when earth is the fastest in its motion around the sun.
d) Motions able to change equatorial coordinates
Unfortunately equatorial coordinates will never keep same data (it would be so easy), They will change depending on 2 main periodic motions called precession and nutation. Precession is because of tide effects on earth rotation. It's a slow motion with a period of 26000 years roughly but with a high range because it can increase until 23.5 degrees roughly (able to change a star's declination until 23.5 degrees ). Nutation (here because of tide too) has a short period of 18.6 years with a weak range compaired to precession : 9 seconds only. So for our programmation, and also to simplify our work it's possible to forget nutation, we lose some accuracy at the end but it will be enough…
We will use those formulas for our work :
–General precession=p=50.290966 0.0222226T” per year, T elapsed time in julian centuries (36525 days) since J2000.0.
-p=m*cos(E)+n*sin(E), E average ecliptic's inclination (23°26’21.448” at J2000.0). m is precession for right ascension, n is precession for declination.
-p=5029.0966 2.22226T-0.000042T^2” per century.
-m=4612.4362 2.79312T-0.000278T^2” per century.
-n=2004.3109-0.85330T-0.000217T^2” per century.
Precession causes on earth a kind of motion “spinned ball with effect” compared with sun. Precession causes a polar star change too; currently it's polaris which indicates north but it will not be always that case.
d) Other things which bring a lot of corrections
We will find atmospheric refraction of course ( Snell-Descartes laws), and Bradley's aberration too, parallax, watcher's height compared to the sea level, the half-diameter's mistake (of our watched object). half-diameter's star (see right height) is forgettable because a star is very little seen from our position on earth; we can forget parallax too because of the very high distance between stars and earth, so parallax is very weak here. we can forget Bradley's aberration too…Fortunately at least we will keep the other things to have some results !
Bradley's aberration : This phenomenon exists because of limited light's speed (roughly 300000 km/s in vacuum) and earth's rotation around the sun (roughly 30 km/s , depends on its elliptic course). It makes in 1 year a slight rotation of all the stars around their “average position” on the celestial sphere; In fact a star is not motionless but turns around a little ring or ellips. Difference between the average position and apparent position is 20.5 seconds… It's not a high difference that's why for our work we don't care about it.
This difference is explained with those draws just above : tan(alpha)=v/c ; but on earth v=30 km/s roughly and c=300000 km/s. So tan(alpha)=0.0001, so alpha=0.00573°=20.6265”.
Little stories :
.It's difficult to see this phenomenon in his mind, so “a little draw is sometimes better to grab every slightest detail of it”, an analogy of situation can be handy. There is Bradley's analogy himself due to his navigations on Tamise, watching a flag onboard : he has noticed that flag wasn't nor in wind direction, , nor in boat's course but in a vectorial addition of them 2. It's the same for a star's position depending on earth 's speed around the sun and light's speed.
. Parallax can be used like an efficient tool to know approximatively distances between earth and planets, stars. For the closest stars parallax is under 1 second…So we will store and forget that it's not really important for us.
Atmospheric refraction's correction it changes depending on atmospheric pressure and temperature ; on wikipedia we find a formula to use when we are on the sea level, with a steadily atmospheric pressure of 1010 mbar with a temperature of 10 °C. This formula, is the same as the formula's master (it's the same but simplified) is written like this :
R=(1.02/tan(Hv 10.3/(Hv 5.11)))/60 ; R is refraction angle in degrees, Hv is the true height in meters.
apparent depression line, due to the watcher's height compared to the sea level is easily findable on a lot of web pages, written like this :
D=(1.77*√(E))/60 , D horizon's decreasing in degrees and E watcher eyes rising in meters.
e) How can i measure a star's height ?
It's the same thing than the height line with the sun. we have to use the vernal's hour angle , with it we will be able to win the star's hour angle ahl :
ahl = ahgy – g + AV (g is estimated point's e longitude, AV is verse ascension).
verse ascension it's the angle between vernal point's angle and star's watched meridian measured like a hour angle, from east to west and on 360 degrees ; it's used on astronomical navigation to point a star.
Right ascension it's the angle beetween vernal point's meridian and star's meridian, measured from west to east and written in hours. it completes 360 degrees of verse ascension. for exemple verse ascension in january 2016 for Betelgeue is 270.985 degrees . so right ascension will be equal to 360-270.985=89.015 degrees. the most of time when we are seeking about a star's coordinates on the web it's given with declination and right ascension, at J2000.0.
Clike the right height with the sun (see right height), formulas to handle for our calculations are those ones :
-sin (hc) = sin (L) * sin (D) + cos (L) * cos (D) * cos (AHL) ; hc is star's height calculated, L is estimated point's e latitude, D is star's declination, and AHL is local hour's angle.
-cos (z^) = (sin (D) – sin (L) * sin (hc)) / (cos (L) * cos (hc)).
– the star's bearing calculated is equal to z^ if 180°<=AHL<360°, and 360-z^ if 0<=AHL<180°.
-Intercept is difference between star's true height seen with sextant, and star's calculated height. Intercept = hv-hc. Intercept is written in miles.
For people who have resist until now in spite of this mess of formulas and unbearable calculations , a little funny story before the next part :
2) Our program
This program works for Dubhé which belongs to the great bear, but it's also able to work for any other star. you have just to replace values t (right ascension of the star at 01/01/16 at 00:00 utc) and u ( star's declination the 01/01/16 at 00:00 utc). this program works from 01/01/2016 at 00h00 utc to the year 2100 which is not a bissextile year (but this program considers it like this, and that is wrong).
1) part 1, change the date and hour in one workable number (line 1/line 36).
Until the line finished by , g (line 36), program will change hour, minutes, seconds, months and years of the watch in a single number called g. to find bissextile years we use a cosinus fonction (a periodic sinus fonction would have been the same); this fonction is equal to 1 every 4 years. When this fonction here called i= 1, we can see that e (day's number depending on the month) changes. Part with “elseif” find the ellapsed day's number in the month since 01/01 of watch's year until the moment of the watch.
int(365.25*f) find year's number between measurement's moment and 01/01/2016. int() is a fonction able to keep the integer part of a number, for exemple int(17.8)=17.
line 35 carves all our parameters in the same unit (the day), They are added up to find g. 2457388.5 (=2451545 5843.5) is the Julian day the 01/01/16 at 00:00.
2) part 2, calculate vernal point's hour angle ( line 37 / line 45).
Line 37, 100+5.5/60 is vernal point's hour angle ahso the 01/01/16 at 00:00 utc.
Lines 37 and 38, with g, we are able to calculate ahso when it has been measured . We know that the vernal point runs in 24 hours more than 360 degrees ; to have its average speed with accuracy we have to divide the distance ran (360 degrees of revolution bitween 2 equinoxes in 1 year) with time made to vernal's point to do it. This duration is the tropical's average year equal to 365.2421898 days. So (360/(365.2421898)+360)/24=15.0410686… In fact the average vernal speed per hour (for 2000. Elle decreases slowly every year).
Currently the average tropical's year decreases about 5,24 seconds per century, -0.0524*16/(24*3600*36525) is calculation for this decreasing time between 01/01/16 and 2016.
-0.0524*(g-2451545-5843.5)/(24*3600*36525) is calculation for this decreasing time between 2016 and the watch's moment.
Between lines 38 and 43, our program gives z (ahso) correct, because before this step z will be quickly higher than 360.
z>DMS changes z writting from decimal to sexagesimal shape (degrees , minutes, seconds).
3) Part 3, (line 46/line 76) calculated watch and intercept.
At the beginning on this step we change sexagesimal data of estimated point and measured angle in decimal data workable for the next calculations. In this program, i could have insert collimation and sextant's excentricity too but i haven't done it (it's a mack 3 model…). Lines 56 and 57 it's right ascension and Dubhé's declination the 01/01/16. Line 58 to 61 we use precession's equations to get Dubhé's celestial coordinates at the measurement's moment. line 62 caculates star's verse ascension at measurement's moment, and…for people who are still reading this article, still able to stand that a little break they deserve it ! (if it could also attract more visitors to the blog…)
So where are we now… well yes ;
program's rest uses formulas seen at the beginning. if fonction uses between lines 70 and line 72 gets true azimut. and now, the test !
3) Test
to manage a star's point we need several things :
-clear horyzon.
-An adjusted clock on UTC.
-A TI89, a programmable calculator or ephemeris.
-A sextant or a tool able to measure an angle.
-An offshore chart.
For me day of the test horyzon wasn't really clear, and i've used an old mack 3's sextant 3 so it wasn't very accurate. besides this it has been done on the english channel, so range of tide could have grown this lack of accuracy.
conditions have let me the allowance to do 3 measures, First on Dubhé, second on Polaris and third on Alkaid (handle's end on the great bear's shape). so finally we stay on 3 close's stars only, but at least it gives a little idea after all. estimated point e (chosen by any chance on the chart) used is 49°50′ North, and 1°5′ East. so e coordinates are positive (negative latitude if it's south, negative longitude if it's west).
Résults:
.For Dubhé, true height is 35°52′ watched at 22h 32 min 52 seconds (local hour, to get it utc we substract 2 hours) the 15/08/16. Program calculates that :
calculated azimut = 325.5° ; intercept = -22.423 miles.
.For Polaris, true height is 48°46′ watched at 22h 35 min 2 seconds. So we get :
calculated azimut = 359 007 ° ; intercept=-47.113 miles.
.For Alkaid, the true height is 44 ° 31′ watched at 22h 28 min 42 seconds. So it gives :
calculated azimut = 299.38 ° ; intercept=49.285 miles.
So when we use the map, we go from our estimated point e, then we draw our intercept ; from e to the star if it's positive, and from e to the same direction but the other sense if intercept is negative. then we can easily draw height lines, perpendicular with intercepts to define our triangular position. So our get point locate us at 80 km southward of Fecamp city (my chart was to small to draw my height lines)… and measurements has been made at Pourville next to Dieppe city… make this again with a clearest horyzon and a better sextant maybe, with more distant stars i guess.
This program works well for clestial coordinates corrections, for vernal's hour angle it's accurate at /- 1° (it makes a mistake's range of 0.3 % roughly). It works with an average tropical's year but it changes every year depending on a lot of parameters , like moon influence, the other planets. vernal's point speed is not regular but changes a little.
To improve it we could include collimation, sextant's excentricity, also doing things to get directly coordinates…have a better accuracy on vernal's point and take care about nutation, Bradley,captain's age…